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I. INTRODUCTION 



o 

Free random variables [l], play an increasingly important role in mathematics, physics, multivariate statistics 
and interdisciplinary research |3l— TlOl] - The cornerstones of this success are the so-called R and S transforms. The R 
transform allows one to infer the spectral properties of the sum of random operators, provided the individual spectral 
measures are known for each of them and they are independent in the noncommutative sense also known as free. 
The S transform plays a similar role for the multiplication of free random operators. These constructions allow for 
fast decomposition of several problems for complicated random operators into simple ingredients. Since free random 
operators have an explicit realization in terms of infinitely large random matrices, the techniques based on the R and 
S transforms provide a powerful tool to solve technically involved problems in random matrix theory in an easy way 
when traditional methods break down. 

Historically, the R transform was devised for hermitian operators and the S transform for unitary ones. The issue of 
the generalization of these constructions to other classes of operators was a subject of intensive research during the last 
two decades. In particular, one of the most challenging problems was the question of the possibility of an extension 
I 1 of the R and S transforms to strictly non-hermitian matrices, which find nowadays vast applications in many fields of 
, research. This problem is also especially interesting as traditional techniques developed for hermitian random matrices 
CNJ 1 generally fail in the non-hermitian case. Some time ago, two of the present authors have extended the additive R 
transform for the non-hermitian ensembles [HKIIj]- Similar constructions were also proposed independently in (ljjLIXij. 
and were soon generalized (l5l [l6| . The question of defining the multiplicative S transform for non-hermitian matrices 
was however open and frequently doubts were expressed whether such a construction is possible at all. On the other 
hand several complicated problems involvi ng p roducts of large matrices have been solved using other methods and 
results were sometimes surprisingly simple [l7H22j |. hinting at the possibility of a hidden mathematical structure. 

In this work we demonstrate that such a structure - the non-hermitian S transform - exists and can be used as 
a powerful algorithm for solving the spectral problems of various products of random matrices. As a byproduct we 
also generalize the ordinary 'hermitian' multiplicative technique to matrix ensembles with vanishing mean which was 
\T2 • never done before. 

In Section[TT]we outline main results of the paper. In particular we give the multiplication law for free non-hermitian 
C$ . matrices. 

In the next two sections, in order to make the paper self-contained, we introduce diagrammatic techniques which 
will be the main tool for deriving the key results of this paper. 

In Section [Til] we very briefly recall the formalism to calculate the eigenvalue densities of large random hermitian 
matrices in the limit of matrix dimensions N —> oo. We recall the connection to planar diagrams and use the 
diagrammatic technique to give a simple proof of the addition law. 

In Section IIVI we repeat the discuss for non-hermitian matrices. We show that the Green's function and the R 
transform are given by 2 x 2 matrices and recall the formalism to handle this case. 
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In Section [VI which is the main section of this paper, we first rederive the multiplication law for hermitian matrices 
using diagrammatic arguments and then we generalize the construction to non-hermitian matrices. We discuss the 
S transform for this case and show that similarly to the nonhermitian versions of the R transform and the Green's 
function it has a form of a 2 x 2 matrix. 

Finally in Section IVII we give examples of application of this law to practical calculations of the eigenvalue density 
for a product of free matrices. We conclude the paper with a short summary. 



II. MAIN RESULTS 



In this section we shortly summarize the main results of this paper. The key quantity of interest in random matrix 
theory is the eigenvalue density, which may be equivalently expressed through the Green's function. The R and 
S transforms satisfy functional relations with the Green's function and hence their knowledge is equivalent (in the 
hermitian case) to the knowledge of the eigenvalue density (or more precisely of its moments) . 

Explicitly, the standard form of the multiplication law of free large hermitian matrices is given in terms of the S 
transform jl| just through an ordinary product 

S AB (z) = S A (z)S B (z) (1) 
The S transform is a complex function of a complex variable and it is related to the R transform as follows 

S ^ - wm (2) 

The two relations given above hold only if matrices A and B are not centered: (trA) ^ and (trB) ^ 0. This means 
the corresponding R transforms may not vanish at the origin of the complex plane Ra(z = 0) ^ and Rb{z = 0) =/= 0. 
If either Ra(0) — or i?s(0) = but not both, the corresponding S transforms do not exist, but one can still save 
the multiplication law [18|]. The prescription [18[ breaks down when both means ( i.e. for A and B) ensemble vanish. 
One of our main new results is that one can still write down a multiplication law in terms of the R transform in that 
case too, using the the following set of equations 

Rab {z)=Ra (w)R b (v) 

v=zRa(w) (3) 
w=zRb(v) 

which involves three complex variables z,w,v. One can eliminate w and v for given Ra and Rb to obtain Rab{z). 
This set is equivalent to the standard equation ([1} when the matrices A and B are not centered but it is also valid 
when either of the two matrices, or even both, are centered, making this a more general formulation. This set of 
equations is quite handy in practical calculations too. One can use it to directly calculate the R transform of the free 
product avoiding the determination of all auxiliary functions and the S transform in particular. Another advantage of 
these equations is that they can be generalized in a natural manner to the case of free multiplication of non-hermitian 
operators and thus they can be used to determine the eigenvalue distribution of products of non-hermitian matrices 
taken from independent random ensembles in the large N limit. 

Before we write down the corresponding set of equations let us first recall that the Green's function for non-hermitian 
matrices are conveniently expressed as two-by-two matrices with complex elements [Til [l2j ■ This will be in detail 
explained in the paper. The R transform in this case is a map of a space of two-by-two complex matrices onto a space 
of two-by-two complex matrices Q — > 1Z(Q). In order to distinguish this situation from the hermitian case ([3]) where 
functions and their arguments were complex numbers we shall use calligraphic letters to denote the corresponding 
two-by-two complex matrices. The law of free multiplication for non-hermitian matrices reads 

n A B{Q)=[n A {QB)] L ■ [n B {Q A )] R 
[Q A f=Q ■ [k a {Gb)] l (4) 

[Qb] L =[Kb{Qa)] R ■ G ■ 

It has almost an identical algebraic structure as © except that now all objects are two-by-two matrices and thus 
the order of multiplications matters. The superscripts R and L outside the square brackets, which were absent in 
([3]), stand for right or left rotations, respectively, of a matrix X in the brackets: [X) L — UXW and [X] R = U'XU . 
The matrix U is a unitary diagonal matrix U = diag(e"^ 4 , e~ Ic ^ 4 ) that depends on the phase <j> of the complex 
number z — \z\e 1 ^ being the argument of the Green function Q = Q(z,z) containing the information on the spectral 
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distribution of complex eigenvalues on the complex plane z. Although this set of equations is more complicated than 
for hermitian matrices (O it also gives a direct, practical way of determining the Green's function Q of the product 
of random matrices A and B. We will illustrate this by an explicit examples towards the end of the paper. We will 
also introduce the S transform for non- hermitian matrices and use it to rewrite the set of equations (U), however we 
think that from the operational point of view equations Q are more convenient. 



III. HERMITIAN MATRICES 



A. Preliminaries 



We are interested in finding the distribution of eigenvalues Ai, in the limit when N (the size of the matrix H) is 
infinite. The average spectral distribution reads 



p(A) = Jta -(EW-A^ (5) 



\i=l 



where Aj are eigenvalues of a random hermitian matrix H and brakcts (. . .) denote averaging over a given ensemble 
of N x N random hermitian matrices generated with the probability 

P(H)<xe- NTtV W. (6) 

For hermitian matrices eigenvalues A^'s lie on the real axis. It is convenient to introduce a complex- valued resolvent 
(Green's function) 

G(z)= hm I /■&_-!_ \ . (7) 



N^oo N \ zt - H 
from which one can reconstruct the spectral density function ([5]) 

p(X) = 1 ^ Urn (G(X-ie)-G(X + ie)). (8) 

using the well-known formula x ± i0+ — P-V.4 T iTrS(X). The symbol 1 will be used throughout the paper to denote 
identity matrices of different size. Here it was an N-by-N identity matrix. The Green's function is a generating 
function for spectral moments /i n = limjv-s.oo jt (TrH n ) 



n=0 

with fj,Q = 1, as follows from the 1/z-expansion of (JTJ). Another fundamental quantity is the "self-energy" S = E(z) 
defined as 

It is related to the Green's function by an independent equation 

X(z) = R(G(z)) , (11) 

where the function 



R(z) = J2*nZ n - 1 (12) 



71 = 1 

is the generating function for planar connected moments n n — limjv^oc ^-((Tri?™)) called free cumulants and denoted 
by double brackets. This function is usually referred to as the R transform. Its form can be deduced from the 
integration measure ((SJ). The difference between the planar connected moments (free cumulants) n n in (|12[) and the 
spectral moments fi n (|9]) will be explained in the next section where a diagrammatic interpretation of these equations 
will be discussed. 
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FIG. 1. (Left) An example of a diagram contributing to the generating function Gij(z). Two end-points should be labeled by 
indices ij. Each horizontal dashed line corresponds to ^5 a t, while a double line represents the expectation value (the propagator) 
{H a bH c d)o = jf^SadSbc- Since all lines in the diagram are proportional to the delta function this equation reduces to a scalar 
equation for G(z). Each horizontal dashed line corresponds after this reduction to ~, each double line to — , each vertex to g n . 
The shown diagram contributes to the seventh moment {TrH 7 } which of order \ in the series expansion ((9| since it has eight 
horizontal lines. The diagram contains seven cubic vertices g\ and one quartic vertex g^ that are generated by the perturbative 
expansion of the residual part of (| 15f) . Each pair of dots on the horizontal line corresponds to a factor H a t> inside the average 
(Tr_ff 7 )o- (Right) The graphical notation for the generating function G(z). It generates diagrams having two end-points which 
include for example the one shown on the left. 



The relation between the generating function for spectral moments G(z) and the generating function for connected 
moments R{z) can be made explicit if one eliminates £ from (fTTJf and (fTTj) . This yields a relation 



which is equivalent to 



G(z) = J^WW) (13) 



G(R(z) + -) = :. (1.4) 



One can use these relations to determine G(z) for given R(z) or vice versa. To give an example, consider the simplest 
case of a random matrix from the Gaussian Unitary Ensemble (GUE). In this case the only non- vanishing cumulant 
is K2- Without loss of generality we can choose K2 = 1, so that R(z) — z. Using (|13p we have G(z) = l/(z — G{z)). 
The last equation can be easily solved for G(z) and the solution can be used to calculate the spectral density 
One recovers the Wigner's semicircle p(X) — — A 2 |23j . 



B. Planar diagrams 

One can calculate (O by Gaussian perturbation theory. One does it by splitting the integration measure (0 into a 
Gaussian part and a residual part 

The Gaussian part Pq(H) is then used to calculate averages (. . .)o while the remaining expression is left inside the 
brackets and is averaged with respect to Pq. The constant Af is an overall normalization. This non-Gaussian part 
is perturbatively expanded in g n , so effectively one has to calculate averages of various powers of H with respect to 
the Gaussian measure. Each term in this expansion has a graphical representation, similar to Feynman diagrams 
known from quantum field theory (sec figure [I}. For example, single horizontal lines represent contributions from the 
factors -t in In the large N limit only planar diagrams contribute to G(z), since all others are suppressed by 
0(1/N) factors (note that each closed line generates a factor N coming from contraction of indices 5u = N). Thus the 
calculation of G(z) amounts to summing all (infinitely many) contributions from planar diagrams with two endpoints 
as shown in figure ([T]). Actually in the most general case one should rather consider a matrix form of the Green's 
function (G = (Gij(z)) where i and j are indices of two end-points i = 1, . . . ,N, j = 1, . . . ,iV (see figure [J) and 
calculate the scalar function (JTJ) afterwards as the normalized trace G(z) = -^Tr<G(,z). Also the self-energy equation 
(|10j) should formally be written in a matrix form. However in our case all generating matrices are proportional to 
Kronecker delta functions z%j = zdij, Gij(z) = G(z)5ij, £.y (z) = S(z)^, (HijHki)o ~ SuSjk so all equations like (fTTJ)) 
and Ijlip reduce to scalar equations for the coefficients multiplying the delta functions. 




•- G 



FIG. 2. Diagrams G(z) can be obtained from one-line- irreducible diagrams (see figure[3]) by joining them one after another. 





FIG. 3. (Right) An example of an one-line-irreducible diagram. (Left) The graphical notation for the generating function £(2) 
of one-line-irreducible diagrams. 



A graphical interpretation of equation (|10[) becomes clear if one rewrites it as an infinite geometric series 



G(z) = - + -£(*)- + -^{z)-^(z)- + ... . 

z z z z z z 



(16) 



which can be seen in figure [3] This figure tells us that all diagrams in G(z) can be constructed by lining up one- 
line-irreducible diagrams one after another. An example of such a one-line-irreducible diagram contributing to £(z) 
is shown in figure |3l Such diagrams are characterized by the property that they cannot be disconnected by cutting 
one line as opposed to diagrams generated by G(z). Indeed, as one can see in figure [5] a diagram in G(z) can be 
disconnected by cutting any horizontal line like that between two consecutive E's. The diagrammatic equation in 
figure [5] can be interpreted as a definition of the generating function S(z) of one- line- irreducible diagrams. 

It turns out that one can write down an independent equation relating S(z) to G(z). One can namely observe that 
any one-line-irreducible diagram can be obtained from diagrams generated by G{z) as shown in figure [4] by adding a 
spider structure making them one-line-irreducible. Each bubble n n of the spider with n double legs corresponds to a 
connected moment (free cumulant) of order n. This equation tells us that 



E(z) = ±((TrH)) + ^((Tri/ 2 »G(z) + l(( Trj ff 3 ))G 2 (z) 



R(G(z)) 



(17) 



The diagrammatic equations in figures [2] and 2] belong to the category of Dyson-Schwinger equations known from 
quantum field theory. They are equivalent to the equations (fTU)) and (fTTj) discussed in the previous section. 




FIG. 4. A one-line-irreducible diagram can be obtained from a one-line-reducible diagram by adding to it a minimal diagram- 
matic structure complying with the measure ()15[) which makes it one-line-irreducible. Such a minimal structure is provided 
by diagrams corresponding to planar connected moments K,k (free cumulants) (see figure [5j) which we indicated by bubbles 
surrounded by double circles in the figure. This double ring around the bubble is chosen to make it similar to double brackets 
used in our notation for connected averages. Diagrams in such a bubble are connected. The difference between diagrams 
corresponding to planar connected moments (cumulants) and spectral moments is explained in figure [5] 



FIG. 5. (Left) An example of a diagram generated by fifth free cumulant K5 = -^((TtH 5 )}. All diagrams in the bubble must 
be connected in contrast to the diagrams generated by spectral moments. (Right) An example of the decomposition of some 
diagrams generated by the fifth spectral moment (15 = jj(TrH°) into two connected moments Some other diagrams in 

/U5 can be decomposed into K1K2K2 or any other combination of cumulants as long as the number of external legs is five. Only 
a small subset of diagrams in ^5 corresponds to those of k$ . 

C. Addition law: R transform 

The R transform [lj is important because it allows one to concisely write down a law of addition of (free) independent 
large matrices. Consider first a factorized measure for two large matrices A, B in the limit ./V — > 00 

P(A,B) = P A (A)P B (B) (18) 

where P A (A) ~ exp -NTrV A (A) and P B (B) ~ exp —NTtVb (B) . Then consider a matrix H = A + B. The law of 
addition tells us how to calculate the spectral density of H for given spectral densities of A and B. 

The idea is based on the observation that connected planar moments (free cumulants) of the sum H = A + B split 
into two independent parts 

±-{{Tr(A + B) k )) = ±{{TrA k )) + ±{{TrB k )) . (19) 

The reason for this separation of connected moments can be easily understood in terms of Feynman diagrams. All 
mixed connected moments ((Tr A a B b A c B d . . .)) disappear just because there is no direct line in a connected diagram 
between a vertex of type A and B since the Ai?-propagator is zero (AijBki)o — 0. The crossed pairs of double lines 
corresponding to A and B vanish in the large N limit, since they represent non-planar contribution. So all external 
lines of a bubble generated by A;-th cumulant correspond either A or to B. In other words free cumulants fulfill a 
simple equation 

KA+B.n = KA,n + ^B,n (20) 

and thus 

Ra+b{z) =R a {z)+R b {z) . (21) 

The argument given above is equivalent to a reasoning based on non-crossing partitions used to prove this law in ■ 
The law of free addition (|2ip is also sufficient to calculate spectral moments of the free sum fi n = -k(TrH n ) — 
■h{Tt(A + B) n ) if one knows the spectral moments of A and B. The recipe follows from the relations (fTB")) and p4|): 

1. Using ([T3)l calculate R A {z) for given G A (z) and Rb{z) for Gb{z). 

2. Construct the R transform R a+ b{z) for the sum using the addition law (|2"T1) . 

3. Calculate G a +b(z) for R a+ b(z) using (fH|) and calculate spectral moments (Tr(A + B) n ) and the spectral 
density of A + B using ©. 

IV. NON-HERMITIAN RANDOM MATRICES 
A. Preliminaries 



We now briefly recall how to calculate the spectral density of non-hermitian random matrices using generalized 
Green's functions (l2| . The crucial difference between the hermitian and non-hermitian case comes from the fact 
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that in non-hermitian random matrix models eigenvalues do not lie on the real axis. In the large N limit they may 
accumulate in two-dimensional domains in the complex plane and the corresponding eigenvalue density 

P(z,z)=l™ c -^£* (2) (*-A^ (22) 

may become a continuous function with an extended support in the complex plane. In particular, in stark contrast 
to the hermitian case, the moments p, n = -k>(TrX n ) no longer determine the eigenvalue density. If one wants to apply 
the Green's function formalism for (|22[) one has to find a representation of the two-dimensional delta function and 
not as in the previous section of one dimensional one ([5]) . A natural candidate is 

<5 (2) (z — Xi) = — lim f . (23) 

With help of this representation one can write 



N 



or 



p(z, z) = lim lim / — -r-s ; tttttt ) (24) 

^ ' ' e^ON^oo^Nf^ (e 2 + \z - X t \ 2 ) 2 / v ' 



, . ld 2 F(z,z) 



where 

JV 

< -u ,\'— x \ .V 



F(z,z) = lim lim ( — V In (\z - XA 2 + e 2 ) ) (26) 
or equivalcntly 

F(z,z) = lim Jim ^Trln ((zl - A)(zl - A f ) + e 2 l)^ . (27) 

One can interpret (|24p as a Poisson equation for electrostatics where p(z, z) is a two-dimensional charge distribution 
and F(z, z) is a electrostatic potential [24T - [26T [ . One can further exploit the electrostatic analogy by introducing the 
corresponding electric field which is equal to the Green's function 

G(z, z) = —— = lim lim / — Tr^ ttt ; ^ — r\ . (28) 

v ' ; dz e^0N->oo\N {zl~X j t)(zl-X) + e 2 t)/ K ' 

up to a coefficient. F is a real function on the complex plane, so it is a scalar field from the point of view of two- 
dimensional electrodynamics while G is a complex function and a vector field, respectively. The Poisson equation can 
be rewritten as a Gauss law in two-dimensions 

p(z,z) = -dgG(z,z) . (29) 

TT 

In the large N limit when the eigenvalues A.; of the random matrix coalesce in a certain region of the complex plane, 
the Green's function G(z, z) is no longer holomorphic. Actually as one can see from the Gauss law (|2l?)) the eigenvalue 
distribution p(z, z) is related to the non-holomorphic behavior of the Green's function. 

Let us make a few general remarks about the way we shall use this electrostatic interpretation. In electrostatics one 
usually applies the Gauss law to determine the electric field for a given charge density. In our problem we proceed 
in the opposite direction. We first calculate the Green's function (electric field) and then we use it to determine the 
eigenvalue density. Secondly, in order to calculate the average ([2"5)l one has to take a double limit. It is important 
to take it in the correct order: first to send N to infinity and only then e to zero, since if one took this limit in the 
opposite order by first setting e = for a finite matrix, then the expression in the brackets in (|28|l would reduce to 
l/NTi(z\ — A) -1 . Finally, whenever we apply generating functions for planar diagrams we can automatically take 
the limit e — > 0, which trivially amounts to setting e = 0, since the large N limit (N — > oo) has already been taken by 
the planar approximation used to write relations between generating functions for planar diagrams. 

Note that the Green's function (]28[) is a complicated object which does not resemble its hermitian counterpart - in 
particular we cannot just apply the geometric series expansion that was crucial for calculations in the hermitian case 
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© ■ We can however use a trick, invented in I2| , which allows us to apply the geometric series expansion but for an 
extended 2N x 2N matrix: 

where we have introduced the block-trace operation 

Trb2 (^ £L xa ;(£* (3i) 

which reduces 2N x 2N matrices to 2 x 2 ones. The elements of Q read explicitly: 

Gii{z,z) = ( — Tr- 

£ii (M) = ( jv Tr 



(zl 


-X^)(zl-X)- 


he 2 l 




-iel 




(zl 


-X)(zl - Xt) - 


|-e 2 l 




-iel 




(zl 


-Xt)Ol-X)- 


he 2 l 




zl~X 




(zl 


-X)(zl-X^)- 


Fe 2 l 



In all these equations we tacitly assume the averages in the right hand side to be calculated in the double limit: first 
N — > oo and then e — > 0. The indices 11, 11 11 and 11 merely reflect positions of blocks in the 2x2 matrix Q. We 
see that the upper-right Gn is equal to the Green's function G(z, z) — Gn(z,z) (|2"5|) . On the other hand, the main 
advantage of using the matrix Q is that it can be calculated using simple geometric series expansion. Indeed, defining 
2N x 2N matrices 



zl iel 
iel zl 



(33) 



and 



(34) 

we can see that the generalized Green's function is given formally by the same definition as the usual Green's function 
G but in the space of doubled dimensions 

Q(z,z) = lim lim -j- ( Tr b2 „ 1 nl ) ■ (35) 

For the sake of the argument we have written now the double limit explicitly. As in the hermitian case, the Green's 
function is completely determined by the knowledge of 'generalized' moments. They are now however matrix-valued 

lim lim — (Tr b2 Z~ Y UZ^U . . . Z' 1 ) (36) 

and are not easily related to the eigenvalue density. As before, we now proceed by applying the diagrammatic 
techniques to determine the non-hermitian Green's function. We begin by writing equations for generating functions 
for planar diagrams. In analogy to (I10p . we introduce the self-energy £ but now as a matrix- valued function 

S(z,zW^f'!j y lli ( Z,2 l) (37) 

As in the hermitian case £ is a generating function for one-line irreducible diagrams. In general it is not diagonal. 
Formally it is related to the Green's function as 

G(z,z) = (Z-E(z,z))- 1 . (38) 

where Z is a diagonal 2x2 matrix 

Z =( 7 0z) W 
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obtained from Z t by taking block trace and setting e = 0. This may be done since the equation (13"5)) is already in the 
limit N — > oo. From here on we will set e = in all equations. 

An explicit solution for the Green's function G(z,z) = Qu(z, z) takes therefore the following form 

G(z, z) = v v ■ (40) 

(z - Sii)(z - En) - En En 



We skipped arguments (z, z) of E's on the right hand side to shorten the notation. The non-diagonal terms in ((30)) 
also contain an interesting information (27j . namely their product is equal to the correlator between left ((£i|) and 
right (|i?i)) eigenvectors of X, introduced originally in [28| 



C(z,z) = -GnGn = £ (^(1^(^)6^ {z - Xi)^ (41) 

Since C(z, z) vanishes outside the eigenvalue support, and for typical nonhcrmitian ensembles is nonzero, the condition 
C{z,z) = often provides a convenient equation for the boundary separating holomorphic and nonholomorphic 
solutions of the spectral problem. Indeed, when off-diagonal terms of E vanish equation (|40[) simplifies to that for 
hermitian matrices G = l/(z — En). 

As in the hermitian case we can write an independent equation relating Q and E - a counterpart of (jlip . The R 
transform however is now a more complicated object since it maps a 2 x 2 matrix Q onto a 2 x 2 matrix E: 

Z(z,z)=1l{g(z,z)) (42) 

or in an explicit notation 



E n (z,z) E n (z,z) J I Tin (G(z, z)) K u {G{z,z)) 



(43) 



In order to complete the analogy to the hermitian case we shall now provide a diagrammatic interpretation of the last 
relation. 



B. Planar Feynman diagrams for non-hermitian matrices 

We shall now discuss the diagrammatic method of calculating eigenvalue densities for non-hermitian random ma- 
trices generated by probability measures of the type P{X) ~ exp (— A^TrV~(A, X*]) in the limit N — > oo, which as 
before corresponds to the limit of planar diagrams. We consider potentials given by sums of terms being alternating 
sequences of powers of X and X^ like XX^X 2 X^ .... Such a potential must be hermitian |V(A", X^)] ^ = V(X, X^) to 
ensure that the expression in the exponent is a real number. The first step of the diagrammatic construction it to split 
the measure into the Gaussian part and the residual one P(X) = Pq(X)P t (X) and use Pq{X) to calculate averages 
(. . .)o which can be represented as Feynman diagrams, exactly as for hermitian matrices (115[) . The Gaussian measure 
Po(V) ~ e - NTlV o( x ) i s constructed from a quadratic potential. The most general form of a quadratic potential being 
a real number is TrVb(A') = aTiXX^ + bTr (X 2 + X^ 2 ) with some real coefficients a, b. The coefficients must be 
appropriately chosen to ensure the potential be positive. The expresion is manifestly positive when expressed in new 
parameters a, r G (—1,1): 

P (X) - exp j-A^y~2 (ttXX^ - r^Tr (XX + X^X^ j (44) 
as one can see for example by writing down the corresponding two-point correlation functions (propagators): 

(x ab xl^ = (xl b X c ^ = —5 ad 5 bc , 

2 

{X ab X cd ) Q = (xl b Xl d ) o = r ■ ^6 ad S bc (45) 

The propagators represent elementary building blocks of Feynman diagrams. As for hermitian matrices the propaga- 
tors are proportional to delta functions, so after taking the block trace we can reduce the problem to a 2 x 2 one with 
propagators corresponding to XX\ X^X, XX 7 X^X^ . The crucial step in inferring the index structure of equations 
relating 2x2 matrices Q and E is to use the correspondence between X «-» Tin and X* f-> Hn which follows from 
equation (|34p . Let us do that. The two-point functions (I45[) reduce to propagators represented by double arcs shown 
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11 11 11 11 11 11 



x xt x x xt xt 

FIG. 6. Propagators for non-hermitian models generated by the Gaussian part of the measure (144[) . The are obtained by 
identification X Wn and X' -o- Hn (|34[). This identification induces the indexing of the endpoints marked as dots in the 
figure (Left) {XX^) Q = (ftiiW n )o = ^ (Middle) (II)„ = {HnHn) = ra 2 ; (Right) (X t X t ) Q = (W„« n ) G = ra 2 . 




FIG. 7. Connected diagrams generated by the fifth order planar cumulant {(XX' XXX'}) (the head of the spider). These 
diagrams contribute a factor Q\\Q\\Q\\Q\\ to Eu = 7?-ii(6). 



in figure [BJ The matrix Z^ 1 (I3£J|) generates lines between 11 vertices which contribute 1/z and lines between 11 which 
contribute 1/z, while there are no lines between mixed vertices. Using these elementary blocks we can draw graphical 
equations as those in figures [2] and [4l The only difference as compared to the hermitian case is that they are written for 
2x2 matrices. Each black dot in the diagrams in these figures is ascribed to an index which may assume two values: 
cither 1 or 1. Each pair of neighboring dots on the horizontal line in figure 2] corresponds to X or or to Hn or Hn 
as follows from the assignment (|34p . As an example consider a spider diagram of order five generated in the expansion 
shown in figure [H Each leg of the spider may be attached to X or X', so on the horizontal line we have a sequence 
of these symbols - for instance XX^XXX\ or equivalently / H.\\K\\H.\\K\yH.\\ ([Ml) . The corresponding diagram 
is shown in figure [7J In a shorthand notation the diagram is determined by a sequence of pairs 11, 11, 11, 11, 11 on 
the horizontal line which begins with 1 and ends with 1 so it contributes to Ejj, since the corresponding diagram is 
one-line irreducible. As one can see from the figure its contribution is proportional to G\iGi\G\\G\i- The indices of 
Q bubbles are enforced by indices of the spider legs - they must match the sequence on the horizontal line. 

All such contributions are captured by a matrix valued function 1Z(Q), in this particular case by its element TZii{G) 
which contains contributions generated by sequences beginning with 1 and ending with 1. Each element of the matrix 
7Z(G) may depend on all elements of the matrix Q so this function maps 2x2 matrices onto 2x2 matrices and in 
general is highly nontrivial (|42l) . The exception is the Gaussian case for which the map is linear. 

For the purpose of this paper let us study Gaussian case in more detail. The most general Gaussian ensemble ((44)) 
leads through (|4"5j) to (see figure [5J) 

Let us now constrain ourselves to the so-called Ginibre-Girko ensemble which corresponds to the case r = and a = 1 
in ([5]), so the matrix S reads 



K(s) = ( si! 3! J = ( sii o J < 47) 

where the off-diagonal contributions are analogous to the relation R(G) = G for the hermitian Gaussian ensemble. 
Solving (|38H47p determines the spectral problem for the Ginibre-Girko ensemble. Inserting (14T1) into (I3"51) we get: 



( Gn Gn \ = 1 ( z Gn \ 

\Gn Gn) \z\i-GnGn \Qii * ) 



(48) 
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The equation for off-diagonal element reads 

\z\ 2 -GnG 



Gil = TZvI^k-F- ■ ( 49 ) 



It has two-solutions: one with Q x \ = and the another one with Q x \ ^ 0. The first one leads to a holomorphic Green's 
function G = Q\\ 

G(z) = - (50) 

z 



while the second one to a non-holomorphic (see the upper diagonal component of equation G = Q\\ 

G{z,z) = z (51) 

which gives the following eigenvalue density 

1 d 1 

p{x,y) = -—Gu{z,z) = - . (52) 

IT OZ 7T 

Both solutions match at the boundary \z\ 2 = 1. So we have recovered a known result that the complex eigenvalues of 
the Ginibre-Girko ensemble are uniformly distributed on the unit disc. 



C. Addition law for non-hermitian matrices 



One can actually use exactly the same arguments as for hermitian matrices to deduce the law of free addition for 
non-hermitian matrices. It has a simple form given in terms of matrix-valued R transforms: 

n A+B {G) = K A {G)+n B {G) (53) 

which follows from the fact that all mixed AB propagators vanish and therefore all mixed connected diagrams 
having a line between A and B vanish too. Since such diagrams represent connected moments (free cumulants), e.g. 
■h ((AB 2 j4J AB'')) = 0, we see that the only non-zero contributions come from connected diagrams (moments) which 
either have all A's or all B's. For applications and more details of this generalized addition law we refer to [12Hl4| . 



V. MULTIPLICATION LAW 



A. Preliminaries 



The S transform plays the same role for matrix multiplication as the R transform for addition. Assume that A and 
B are large independent (free) random matrices given by a product measure (|18p . The multiplication law tells us how 
to calculate spectral moments (Tr(Ai?) n ) of the product H = AB provided we know the spectral moments of A 
and B or equivalently that we know the corresponding Green's functions Ga(z) and G B [z). The multiplication law, 
expressed in terms of the S transform, reads [l| 

Sa-b(z) = S A (z)S B (z) (54) 

and the S transform is defined by 

fif( z ) = i±£ x ( z ), where X {zG{z) - 1) = ^ . (55) 

The algorithm for " multiplication" is similar to that for " addition" : 

(i) Calculate S A (z) and S B (z) using (|5"5"|) . 

(ii) Use the multiplication law (|54"1) . 

(iii) Use again ([5"5j) to derive G AB (z) for the product of AB. 

Let us first derive some useful relations between the R and S transforms. Changing variables z = yG(y) — 1 in (|55[) 
we get 



S(yG(y) - 1) = _ X 

y G(y) 



(56) 
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Using (|10p we can rewrite the last equation as 



S (G(y)Z(y)) = ^ (57) 



Setting T,(z) = R(G(z)) and taking the reciprocals of both sides we arrive at 

1 



R(G(y)) . (58) 



S (G(y)R(G(y))) 

Changing variables once again to z = G{y) we obtain the equation 

which gives an explicit relation between the R and S transforms. The S transform can be defined only if the R 
trasform does not vanish at the origin: R(0) ^ 0. This corresponds to random matrices with a non-vanishing first 
moment (cumulant) -^(Trif) = -^((TrH)) ^ 0. Otherwise the S transform cannot be defined as a power series and all 
the manipulations presented above break down. The last equation can be inverted. Let us introduce a new variable 
y = zR(z). Now d^J) reads 

S(y) = / N = —r^ r = — ^ r (60) 



R 



where z can be recursively eliminated by repeating the substitution z = ad infinitum. This leads to a function 

which is nested infinitely many times forming a sort of continued fraction. The last equation can alternatively be 
written as 

s ^ z) = (61) 

which is an inverse formula to (|59p . The two equations can written in a symmetric way as mutually inverse maps 

z — yS(y) and y = zR{z) . (62) 

As an example we consider a shifted Gaussian random matrix which has only two first non- vanishing cumulants. For 
the standardized choice k% = Ki = 1 the R transform reads R(z) = X + z. Using (|6"Tj) we obtain 

S ^ - T+zS(zj (63) 



and hence S(z) = z^+VHM 



2z 



B. Diagrammatic derivation of the multiplication law 

We are now ready to diagrammatically derive the S transform and the corresponding multiplication law. The 
argument given below will turn out to be crucial for the generalization to non-hermitian matrices. The initial point 
of the construction is to consider a 2N x 2N block matrix T-L and its even powertQ 

The upper-left corner of H 2k involves solely the powers of AB, which we are interested in. In order to have an 
access to the traces of individual blocks in the matrix we again apply the block trace operation defined before. The 
upper-left corner of the reduced matrix Trb2% 2fc is equal Tr(AB) k while of Trb2"H 2fc+1 is equal zero. So now the idea 
is to reformulate the problem of calculating the Green's function for the product 



1 This should not be confused with the 2N X 2./V block matrix constructed for the nonhermitian Random Matrix Ensembles in section IV. 
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FIG. 8. The spider diagrams correspond to free cumulants generated by the matrix A = H12 and therefore the double dots 
on the horizontal line are indexed by 12. So on the horizontal line we have alternating indices 1212 ... 12 and this enforces all 
the (J-bubbles to have indices 21 as one can see in the figure. Therefore there are only Q21 bubbles in the diagram and the 
corresponding equation is £12 (w) = Ra(G2i)- The analogous equation for B-cumulants is £21(11;) = Rb(Qi2)- Similarly, one 
can see that £n(u>) = £22 (w) = since one of the double dots on the horizontal line would need to have identical indices, for 
which as we know from l|64p the double line is equal zero. 



as a problem of calculating the upper-left corner of the Green's function Q(w) for the matrx W: 

ri \ ( Qu{w) Qi2{w) \ 1 / 1 \ ( , 

Giw) = { o*U g 22 (w) ) = n \ h2 ~wt^Ti) ■ (66) 

where w is a complex number and 1 is a unity matrix of dimensions 27V x 27V. One can easily check that 

Gab(z = w) = (67) 



since only every second (even) power of Q(w) contributes to the power expansion of Qxi(w), which is thus a power 
expansion inz=w 2 . 

The next step is to define self-energy £(u>). It is a 2 x 2 matrix 

_ / Sn(w) T, 12 (w) 

that is related to the Green's function as 



EH= " ~ (68) 

v ' 1 2j 2 l(w) T>22\W) 1 



G(w) = (wl - £H) 1 (69) 

in analogy to (|TD|) . All the matrices in the last equation are of 2 x 2 dimensions. This is the first Dyson- Schwinger 
equation. To write down the second one - a counterpart of (1111) it is convenient to use its diagrammatic representation 
as that in figure SJ Instead of a scalar equation (fTTj) we will have a matrix equation for 2x2 matrices Q and S (|6"6")l 
and (|68[) . Since we have now have 2x2 matrices it is crucial to work out the index structure of the corresponding 
equation. This structure stems from the correspondence A -H> H12 and B «-» H21 that follows from the position of the 
blocks in H (|64[> . Note the difference to the case discussed in the previous section were we had diagonal blocks (l34l) . 

The only non-vanishing cumulants are ^((Tr"H™ 2 )) = -^{{TrA n )) = n A>n or -^{{TrU^)) = j r (( r TrB n )) = n B , n 
while all mixed ones vanish as we discussed in the previous section. Due to this, the index structure of non-vanishing 
one-line- irreducible diagrams is restricted to that shown in figure [5] and its counterpart obtained by exchanging lf>2 
and A O B. The diagrammatic equations discussed in figure [8] can be summarized as 



Insering this into (|6"9")l yields 



- ' = I »(,)) RAle l llw)] I 



' Gn(w) Su(w) \ _ f w -Ra{Qh{w)) 

v 02l(l") Q22{w) ) \ -R B (Gl2{w)) W 



(71) 



which gives a direct relation between the Green's function Q(w) and the R transform. 

We shall now rewrite this equation in a way which explicitly exhibits multiplicative structure. Note that in the 
following manipulations we do not need to assume anything about the first moment i.e. whether the ensemble is 
centered or not. Inverting the matrix on the right hand side we obtain 



Gi2(w) = —— Ra (021 («?)) , G 2 i(w) = -^—Rb 
Det Det 



(Si* W) (72) 
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Gix(w)=g 22 (w) = ^. (73) 
where Det is the determinant of the matrix, wl — E(u>), on the right hand side of (|71l) : 

Det = w 2 - R A (021 (w)) Rb {Qi2{w)) . (74) 
Inserting two last equations to (167[) we obtain 



r / n _ OiiH J_ = i (7t ., 

AB[ ! w Det z-R A (g 21 (w))R B (G 12 (w)) [ ' 

where z — w 2 . Comparing the denominator in this equation to that of the standard equation G AB (z) = l/(z — 
Rab(G A b(z)) O we get 

RAB(G AB (z)) = R A (g 21 (w))R B (g 12 (w)) . (76) 

At this stage we see the first hint of a multiplicative structure emergence. In order to complete this equation we also 
need (|72[) . Let us set g = G AB (z), g A = g\ 2 {w) and g B = g 2 i(w) to simplify arguments in the R transforms in the 
last equation. Using this substitution we can write (1761) and (I72p in a compact form as a closed set of equations for 
the R transform of the product 

R AB {g) = Ra (.9-B ) Rb (gA ) (77) 

and 

g A = gRAigs) , gB = gRs(gA) • (78) 

which is equivalent to (J3J) announced at the beginning of the paper. This is the multiplication law formulated in terms 
of the R transform. Its main advantage in comparision to the S transform is that it can be applied even to centered 
ensembles (i.e. having vanishing mean) including the case when both are centered (see the examples in sections VIA 
andVIB). 

The difference with respect to the conventional multiplication law S AB (z) = S A (z)S B (z) is that the individual 
factors appearing in (|77p are not expressed uniquely in terms of the properties of a single random matrix ensemble 
e.g. the factor R A (-) is evaluated on g B which is related to the ensemble B. However it is straightforward to obtain 
from (|77|) - ([75|) the conventional multiplication law as we shall illustrate below. 

Let us introduce a new variable y — gR AB (g). We can now express g B - the argument of R A purely in terms of 
the properties of ensemble A: 

u , , Rab{9) y y ,„ n , 

g B = gR B (g A ) = g—— — - = — — - = ^ (79) 

RaK9b) Ra (gB ) Ra 

Now each of the factors in (|77|) depends on a single ensemble. We may do the same for the left hand side, which 
becomes (|60)l 

Putting these formulas together, we can finally write (|TT|) using only the variable y [U 

/: ''u ; | j 1 r|=i^| T 1 r life ( 7 1 rl (^i 



Rab(t^)J \Ra(t^)J {Rb^ 



which amounts to the standard formulation for the multiplication law [l| as follows from (l60l) 

S AB (y) = S A (y)S B (y) . (82) 

The necessity of assuming noncentered distributions comes from the fact that the implicit continued fractions appear- 
ing in (|60l) make sense only for R(z) ~ const + O (z) with nonzero constant term [29(. 



15 



C. Multiplication law for non-hermitian matrices 



In order to derive the multiplication law for non-hermitian matrices we combine the two formalisms outlined in 
previous sections. First we define a 2N x 2N matrix T> in analogy to (|64l) 



V 



A 
B 



(83) 



2Nx2N 



and then duplicate it using (|34|) to obtain an extended Green's function for non-hermitian matrices. This technique has 
been introduced in [2l[ and used for specific ensembles 19, 21]. In this paper we will use it to obtain a multiplication 
law for arbitrary (free) nonhermitian matrices^. 

This procedure leads to a four- fold matricial structure (" double doubling" ) where the primary object is a AN x AN 
matrix 



H = 



and the corresponding Green's function 



G(w,w) 



V 

pt 



r 


.4 







B 

















Q 




V° 





At 


J 



(84) 
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wl 










f ° 


A 





^ 







wl 







B 



















wl 













st 













wl ) 







At 


J 





(85) 



Using the block-trace operation tr (,4 we reduce the problem to calculations for 4x4 matrices 



G(W) 



Gn 


G12 


Gn 


G 


G21 


G22 


G21 


G 


Q11 


G12 


Gn 


G 


G21 


G22 


G21 


G 



■^tr b 49iw, w) 



(86) 



4x4 



where VV = diag(«;, w, w, w). The labeling of the matrix elements follows the convention adopted in the previous 
sections. Similarly, we define a self-energy S(W) as a 4 by 4 matrix: 



G{W) = (W-E(W)) 
which is related to a 4 by 4 matrix representing the generalized R transform: 

E(W) = K(g(W)) . 



(87) 



(88) 



The elements of E and 71 are indexed in the same way as the elements of G Q86p. 

We exploit these 4x4 matrices as auxiliary objects to derive relations between 2x2 Green's functions Ga{Z), 
Gb(Z) and Gm{Z) for A, B and the product M = AB. The naming convention for elements of 2 x 2 matrices 



Ga(Z) = ( G r ^ 11 Q r ^ n 



(89) 



is a bit inconvenient since it requires three subscripts for each element. To avoid multiple subscripts like (A)ll we 
introduce a shorthand notation substituting multiples indices like (A) 11 by AA etc. In this new notation a double 
subscript identifies both the matrix for which the generating function is calculated and the position of the element. 
Using this convention we have 



Ga{Z) = ( G Q AA 

\ y AA VAA 



(90) 



2 To remind the reader, 'free' means essentially that the probability distributions of the two ensembles are independent and that we take 
the N — Y 00 limit. 
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and similarly for two remaining generating functions 

TaA £ 



T A {Z) 



.4.4 



£.4.4 



K A {9) 



( TZ-AA Ti-AA 
V T^AA T^AA 



(91) 



We use the same convention for all matrices, including B and M. For brevity we skipped the arguments of the matrix 
elements on the right hand side of the equations above. We tacitly assumed that they are identical as on the left hand 
side. We will frequently use this shorthand notation below. 

To summarize the notation, TZm denotes a 2 x 2 matrix of the R transform for M while TZmm ~ its upper left 
element, etc. For 4x4 matrices like Q(W), £(W) and TZ(Q) we instead use the indexing as in ([86")) which uniquely 
identifies the positions of elements in such 4x4 matrices. The link between the two conventions emerges from the 
equation 



H = 





H12 


Hn 


u 


H21 


H22 


H 2 i 


H 


Hii 


H12 


Hn 


H 


H21 


H22 


%2i 


H 



B 








.4+ 






fit 





(92) 



that allows us to identify A -f-s- H12, A^ <H- T-L^l and B H21, B* ■<-» H.^. We use this identification to rewrite the 
equation ([88} in terms of 2 x 2 matrices. We begin by noting that even powers H 2k of H (|9"2)l generate powers M k of 
the product M = AB in the upper left corner of the block matrices 

\k 
{BA) k 





n 2k — 



I {ABf 









(AB)t fe 








{BA)V° 



(93) 



These moments are generated by the element Qn(W) of the 4x4 Green's function Q(W) (|86|) or alternatively by the 
element Qmm{Z) of the 2x2 Green's function Qm(Z), so we have 

011 (W) 



Gab(z,z) — Qmm{Z) = 



w 



(94) 



where Z = diag(z, z), W — diag(w, w, w, w) and z 



w 



analogously to (|67[) . This equation, allows us to determine 



Green's function Gab(z, z) and additionally provides a link between Qm and Qa and Qb since elements of the 4x4 
Green's function Q can be explicitly expressed in terms of Qa and Qb, as we will see below using planar Feynman 
diagrams. 

First we recall that all mixed connected diagrams vanish since AB propagators are equal zero. The last statement 
means that there are no direct lines in the diagram connecting A and B vertices. All non- vanishing connected 
diagrams are either of ^4-type like ((■irtrAA'AA...)) or B-type like {{jjtxBB*BB . . .)). They are generated by 
alternating sequences either of A and A^ or of B and B^ but not mixed ones. In other words there are only A-spider 
or _B-spider diagrams. In the T-L notation the first type is generated by sequences of T-L\2 and H^l while the second type 
of H21 and as follows from the correspondance (|9"2j) . We show in figurelHlan example of a diagram contributing to 
the left hand side of equation (J88J. More generally, diagrams with a A-spider have on the horizontal line alternating 
sequences like 'Hi2'H2i'H2i'Hi2'Hi2 . . . which are sandwiched by <J 2 2 , Q12, Glii G21, ■ ■ ■ which match the index sequence. 
The left most index in the sequence of Ws may be equal 1 or 2 and the right most 2 or 1 so the corresponding diagrams 
contribute to S12, Ejj, £22 or £21 • Diagrams with a B-spider have sequences like H21H12 • • -^21 etc, whose left most 
index is either 2 or 1 and the right most index is either 1 or 2, so the corresponding diagrams contribute to £21, £221 
Tiii or £12- AH others E's must be equal zero 



£11 — £22 — £11 — £22 — 
£12 — £21 = £12 — £21 = 



(95) 



since there are no mixed Ai?-spiders. Coming back to the equations for £12, £ij, £22, £21 generated by the A-spider 
we notice that the indices of the Q bubbles which enter the sandwich between the spider legs have complementary 
indices Q22, Q12, Qii, Q21 as compared to those of £'s. The same holds for equations for indices of £?'s and E's 
generated by the B-spider. Moreover, if we compare indices of E's for A spiders to Q's for B spiders we see they 
are identical, and the same holds for E's for B spiders and C/'s for A spiders. All these observations can be concisely 
summarized by the following equation 



(96) 



/ 





Taa 


^AA 







Tbb 








^BB 




^BB 








^BB 


V 





^AA 


£yL4 
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FIG. 9. An example of A-spider connected diagram. Such diagrams are generated by sequences of A and A* which due to 
the correspondence (|92p A H12, A* o H21 generate sequences of pairs 12 and 21. In this example we have a sequence 
12,21,21,12,12 which begins with the index 1 and ends with 2. It contributes to £12 = 7t\2(G), a product of G22, G12, G11, G21 
which can be read off from the picture by matching the indices on the horizontal line. 



The matrix has eight zeros which correspond to (|95l) . The remaining eight elements can be grouped in two groups of 
four elements each of which can be mapped into a 2 x 2 matrix. More precisely, the matrix E of dimensions 4 x 4 is 
expressed in terms of 2 x 2 generating functions 1Z and Q for A and B: 

* A - { *aa Zaa)-{ Kaa(Gb) K aa (Gb) ) Ua{Gb) (97) 



and 



Ebb E bS \ _ / TZbb{Ga) T^bb(G- 



V L ss ^bb J \i<-bb{Ga) R-bbKGa) J 
The argument of TZa in E^ = TZ-a (Gb) is Gb while the argument of 1Z B in Eb = TZb{Ga) is Ga as argued above where 

GA=(t t) ■ GB=(t t) • (99) 



G22 G21 J ' \ Gn G12 

So far we have used diagrammatic properties of the equation ([55]) . Now we can also exploit the second equation (|87[) . 
Inverting the matrix on the left hand side of this equation for the particular form (|96p we can find elements of Q as 
functions of E's. In particular the equation for Q\\ is 

r ww 2 - wE AA E BB - w^AA^BB , lnfV , 
Gl1 = det(W-E) (100) 



where 



dct(W — E) = w 2 w 2 — w Y, AA Y, BB ~ w2 ^aa^bb 

+ V^AA^AA ~ ^AaEaa)(^'BB^'BB ~ ^BB^Bb) 

- ww(T, AA Y, BB + Y, Aa Y, bb ) (101) 



Now we can use to compare Gn/w that follows from (|100j) to 6mm 

" det(Z-S M ) (102) 

where 

= (2 — £mm)(z — E M - M ) — S MM E MM (103) 

and z = if 2 . From this comparison we can deduce relations between Em and E^ and Eb. The numerators of 
expressions for Gn/w and of Gmm are equal if 

^MM — —^AA^BB + ^AA^BB (104) 
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and the denominators (|101[) . (|103[) if 



(w 2 - Emm )( w 2 - E MM ) _ ^mm^mm 

= w 2 w 2 - iu 2 S^Sb5 - w 2 Y, AA 'E BB 

+ (^A^AA ~ ^•AA^Aa)(^BB^'BB - ^BB^Bb) 

- ww(E a a£ Sb + Y, AA Y, BB ) . 

One can check that the two equations are simultaneously fulfilled if 

w 

Emm = ^AA^BB + —^AA^BB 
W 

— ^AA^BB + y ~^pAA^BB 

/w Fw 
— ^AA^BB + y—^AA^BB (106) 
W \ W 

W 

Remarkably, these equalities can be written in a factorizable matrix form as 

V s ™ £ ™/ \\/^aa z AA j \^^ BB e bb j 

= E^E* . (107) 
In order to simplify the notation it is convenient to introduce a 2 x 2 unitary diagonal matrix U 

,1/4 



U = I W , "i/4 I = | „ Jk I (108) 




\ _ / e +l * 

(f) 1/4 ; v ° 

where the angle ip is the phase of w: w = \w\e 1 ^ . Note that w is related to the original variable z as z = u; 2 , so Arg 
z = 2 Arg to. Using this matrix we can associate with any matrix X two similar matrices X L and X R obtained by 
" left and right [/-rotations" of the matrix in question 

X L = [X} L = UXU^ , X R = [X] R = U^XU , (109) 

In particular 

Ei = [E A ] L = UHjJjt , E^ = [E B ] fl = U^ B U . (110) 
The operations [. . .] and [. . .] R obey simple rules like for instance 

[XY] L = [X] L [Y] L = X L Y L , [X- 1 ] 1 = ([A] L ) _1 , X = [[X] L ] R . (Ill) 

which we will frequently use below. 

Now we come to the main result of the paper. Recalling that E^ = 1Z a (Gb) and E^ = 1Z a (Gb) we have (|107[) 

n M {Q M ) = [n A {g B )] L ■ [n B {Q A )] R . (112) 

This equation is a cornerstone of the matrix multiplication for non-hermitian matrices. Let us note the similarity with 
the corresponding equation for the hermitian case ([77]) albeit with two key differences. Firstly, the objects appearing 



in (|112[) are generically noncommuting 2x2 matrices and hence the ordering is crucial. Secondly, the left- and right- 
[/-rotations have no analogue in the scalar hermitian case. 

In fact, to arrive to this point we have only taken advantage of the equations for the element Gn of the 4x4 Green's 
function. Inverting the matrix on the right hand side of (|87[) for E given by (|96p we can relate remaining elements of 
Q to the elements of 2 x 2 E's and 7?.'s. In particular we can write equations for elements Q\2 1 Gil, Q22 an d Q21 which, 
as we know (|99[) form a 2 x 2 matrix corresponding to the Green's function Q A and similarly for C/2i> £/22> ^Ii an d Q\2 
corresponding to Q B . This allows us to express Q A and Q B in terms of E^ and T, B . After some straightforward but 
tedious algebra we arrive at remarkably simple equations which again are analogs of the hermitian equations (|75|) but 
with specific ordering and appropriate [/-rotations 

r riL 



Qm ■ [Tl A {Q B )\ L 



Qb= 



[n B {g A )f-g M . (m) 
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The set of equations (|112p and (II 13[) gives the multiplication law for non-hermitian matrices and constitutes one of 
main results of this work, as mentioned at the beginning of the paper J?]) . 

These equations are in one-to one-correspondence to (177|) and (|78|) except that now instead of complex numbers 
9m, 9a arid gs we have 2x2 matrices Qm, Qa and Qb and the additional [/-rotations. The logic of the method to 
calculate the Green's function for the product M = AB is the same as for hermitian matrices, that is for given Qa 
and Qb one determines the matricial R transforms TZa and TZb and then applies (|112p and (|113D to derive the Green's 
function for M. We will present examples in the next section. Before doing that let us show how these equations can 
be reformulated in terms of a nonhermitian generalization of the S transform. 



D. S transform for non-hermitian matrices 



It is natural to anticipate that the S transform for non-hermitian matrices has a form of a 2 x 2 matrix. It will 
however appear in two different "left" and "right" versions, since 2x2 matrices do not commute in general. To 
demonstrate this, we repeat the arguments which have guided us from ([77)1 and (|78D to (1521 . but now we adapt the 
reasoning to the case of 2 x 2 matrix- valued transforms. 

The first step is to eliminate Qb and Qa from the right hand side of (11 1 2[) and substitute them by Qm in order to 
have the same argument on both sides of the equation. To make the following equations slightly more readable we 
shall skip the subscript M of Qm 

than X^ 1 to avoid too many superscripts. Using (|113l) and (II 1 2[) we have 

R 



writing Q = Qm and we will denote the inverse matrix of a matrix X as -y rather 



Qb = [n%(Q A )Q\ 



1 



-n M {G)Q 



(114) 



This is an equation for Qb but Qb is also present on the right hand side. We can however eliminate Qb by replacing it 
recursively by the right hand side and repeating this infinitely many times. In this way we obtain a nested expression 
(denoted below by dots) 

i R 



Qb = 



1 



n M {Q)Q 



(115) 



that depends on Q and not on Qb 
function of Q: 



Thus we can write the first factor, TZa(Qb), on the right hand side of (|112[) as a 



n L A (Q B ) = n\ 



n M {Q)Q 



i 



SS^> (Km(Q)Q) 



where SS^ is a left S transform defined as 



1_ 

W 7 ^ 



(116) 



(117) 



Let us make two further remarks concerning the notation. In the last equation we skipped the subscript A of SS and 
1Z since the relation is valid for any matrix. The superscript (L) of SS is used on purpose in parentheses to distinguish 
it from L and to emphasize that the left S transform is not a left rotation of the S transform SS^ ^ [SS] L = USU^ 
in contrast to the notation 1Z L = [JZ] L . The function S^ is just defined by the equation above. This equation is 
equivalent to 

ss< L \y) = — 7 ^ (us) 



and 



n L (y) 



n L ([ss( L \y)y] 



ssw \[n{y)y] 1 



(119) 



in analogy to the hermitian case discussed in section IV Al Now we can repeat all the steps for the second factor on 
the right hand side of (|112j) . The result can be written using a right S transform, which is given by two equivalent, 
reciprocal, equations analogous to those of the left S transform above: 

ss( R \y) = — — (120) 



n R ([yss( R )(y)Y 
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K R 0>) = ——7 • ( 121 ) 

ssw {[yn(y)] ' 



Using the left and right S transforms we can write (|1 12[) in a concise form 



-i— = SS^ (G1Im(G)) ■ SS™ (Hm(G)G) (122) 
i<-M{y) 

that depends on G on both sides. This is an equation for the R transform TZm(G) which in turn determines the 
generalized Green's function giving the eigenvalue density. 

Let us rewrite now the left-hand side using either equation (|119j) or (|121[) 



([Tl(G)Gr)\ = [SS { ^> {[GU(G)] U )\ = SS)? 1 (Gn M (G)) ■ SS% ! (n u (G)G) (123) 

which now (almost) takes the form of a multiplication law for SS transforms with the only subtlety being the 
noncommutativity of the arguments. 

In the special case when G and 1Z(G) commute 

[G,n M (G)}=0 (124) 

we get a direct analogue of the hermitian multiplication law for S transforms since all functions are evaluated on the 
same argument y = GTt{G) = 1Z(G)G- In this case it would make sense to introduce yet another S transform: 

n&) = smm ■ (125) 

which does not involve any left or right [/-rotation. It is easy to see that in this case the equation (|122[) can be 
rewritten as 

ss(y) = ssjf } (y) ■ ssf ] (y) . (126) 

One should note that the 2x2 formalism, which has been developed here for non- hermitian random matrix ensembles, 
contains also the standard hermitian case. For hermitian matrices, namely, the Green's functions and the R transforms 
reduce to diagonal matrices 

<127) 

Moreover G = G L = G R , Tl = Tl L = Tl R because the matrix U that defines the left and right rotations (|109[) is 
diagonal too and the product of the diagonal elements gives one. It follows also that SS = SS<- R > = SS<- L '> and that 
the S transform is diagonal SS(Z) = diag(S(z), S(z)) too. Therefore in this case (|122j) takes a diagonal from 



L 



( R ) rr.-v rr.w ec( £ ) 



( S M {z) \_( S A (z)S B (z) 

I S M {z) \ S A (z)S B (z) 



(128) 



that is equivalent to 



VI. EXAMPLES 



In this section we will illustrate our methods by presenting three examples. We will start from two examples which 
cannot be treated, even in the hermitian case, by the conventional S transform treatment as both of the random 
matrix factors of the product are centered. Finally we treat a more complicated example of obtaining a nontrivial 
two-dimensional eigenvalue distribution for a product of two simple factors. 
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A. Product of two free Ginibre-Girko matrices 



Let us first consider the product M — AB of two identically distributed free Ginibre-Girko matrices: A and B. 
Throughout this section we will parametrize matrix elements of the 2x2 Green's functions (|3T)]) with two complex 
functions a = a(z, z) and b = b(z, z) 

(129) 



fan 


Gxx^ 


i f a ib\ 


Uii 


Gn / 


' = {ib a J 



The R transform for a Ginibre-Girko matrix reads (jT7|) . 

a ib \ \ ( ib 



n W= n ((il a l)) = (ib ) (]8 "> 

and its left and right versions 



respectively. We recall that w is related to z as z — w 2 . Let us now apply the multiplication law for M = AB where 
A and B are Ginibre-Girko matrices with unit variance. Using (|112[) we have 

*— i /j., ---II ± ^^r-i^_4 t j (i3 2 , 





ib A 




Since both A and B are identically distributed they have identical Green's function, we can thus reduce the problem 
by introducing a single function b — b A = 

We can now use the two remaining equations of the multiplication law (|113|) which can be conveniently written as 

Qii \Qa] r = [n A {g B )] L 

[Gb] L Gm~ = [Hb(Ga)] r (134) 

In case of identically distributed A and B one of the two equations is redundant and thus it is sufficient to use only 
one of them, for instance the first one. We first eliminate Qm from this equation by using the relation = Z — TZm 
with TZm given by (|133[) : 

(Z - TZ M ) [Ga] R = [n A (G B )} L ■ (135) 





This is an explicit equation for a and b 

(136) 



™ 2 + m 2 



It can be easily solved. It has two solutions: a trivial and a = b = and a non-trivial one a = 0, \b\ 2 = 1 — ww. 
The latter one is equivalent to a — and \b\ 2 = 1 — \/z§ when expressed in the variable z — w 2 . This solution holds 
inside the unit circle: zz < 1 on the z complex plane while the trivial one outside. The boundary of the eigenvalue 
distribution in the z plane is given by the condition 6 = for the non-trivial solution which leads to the unit circle. 
Inserting these solutions to (|133|l and calculating Qm we find 

Gm{z,z)=(S Q \ , for |*| <1 (137) 



Gm(z,z)=( Z q e °A , for |*| >1 (138) 

from which we obtain a rotationally symmetric eigenvalue density for \z\ < 1: 

p(x,y) = -§-_G(z,z) = ^-l- (139) 
7r oz Zir \z\ 

inside the unit circle and p(x,y) = outside. We remind the reader that G(z,z) is equal to the upper left element 
Gmm of Qm ■ 
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B. Product of two free GUE matrices 



We would like to discuss a simple but very interesting case of the product M = AB of two matrices from the 
Gaussian Unitary Ensembles. Both A and B are hermitian but their product is not. Since both matrices have a 
vanishing mean the traditional use of S transform leads to contradiction, as shown in However, our algorithm 
works in this case without any problems. 

Before we apply the full non-hermitian version of the multiplication law let us check what happens if one applies 
its hermitian version given by equations (1771) and (|78[). One can do this since A and B are hermitian. However the 
result for Gm(z) can be interpreted only as a moments' generating function but not as a full Green's function. In 
particular one cannot use it to reconstruct the eigenvalue density ([5]) since the eigenvalues are not constrained to the 
real axis. 

For a standardised GUE matrix we have R(z) = z and thus the multiplication law ([77)1 and (|78D simplifies to 

Rm{9) = 9b9a, 9A = 99b, 9b = 99a (140) 
The two latter relations yield an equation gA = g 2 9A- Its solution is gA = giving Rni(g) = and hence 

G M (z) = ~ (141) 

The moments mt are given by coefficients at l/z k+1 of the 1/z-expansion of Gm{z). We see that all they vanish 
except the trivial one mo = ^(trM°) = 1. Of course it does not mean that all eigenvalues of M vanish. In order to 
determine the eigenvalue density of M one has to apply the full multiplication law in the domain of non-hermitian 
matrices (|112l) and (II 13|) . The calculation goes along the same lines as in the previous example except that now 
instead of f| 130[) the R transform is 

k <«= k ((^))=(^)= s ,142) 

as follows from (I46p for r = 1. It is easy to see that the solution is exactly the same as in the previous example since 
for a = (which was a solution) all equations reduce to those for the previous case. This result is in agreement with 
the recent works 

[HIM HI- Actually one can see that the same holds for any elliptic ensemble with 

since again for a = the equations are identical as before. Again this is in agreement with [li| where it was shown 
that even for A and B being from different elliptic ensembles (ta =/= tb or ta = tb) one obtains the same circular law 

(fim 



C. Pascal limacpn 



We shall calculate now the eigenvalue distribution of the product of two shifted Ginibrc-Girko matrices M = AB = 
(1 + Xa)(1 + Xb) where Xa and Xb are free Ginibre-Girko complex matrices. The main difference to the cases 
discussed before is that the multiplied matrices A and B are not centered: ^trA = 1 and -^tri? = 1, so their first 
moments (cumulant) are not zero: 

SI t A ))-U1 A ) (144) 

Since A and B are identically distributed we set b = bA = bs as in the previous examples. Using (|112j) we have 




(145) 
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Inserting this to (I135[) we obtain an explicit equation 
-ib(*f^ + ,/%) w 2 - 

\\Jw V w I 

which reduces to two equations for a and \b\ 2 : 





(146) 



a(w; 2 -l + |&| 2 -) + |6| 2 (l + -) = l 

w w 

-a(l + -) + -(w 2 -l) + |6| 2 = l (147) 
w w 

This set of equations has a trivial solution: 6 = and a — \/(w 2 — 1) and a non-trivial one that can be found by 
eliminating a from the last set of equations. This gives an equation for C = \b\ 2 (|4"Tj) : 

C 2 + C(l + 2|w| 2 ) + \w\ 4 - \w\ 2 - w 2 - w 2 = (148) 

The border line between the two solutions can be found by setting C = in the last equation [2l[ : 

w 2 w 2 — ww — w 2 + w 2 (149) 

It represents a curve on the z-plane called Pascal's limagon after Etienne Pascal (1588-1651) - the father of Blaise 
Pascal. It has a more familiar form in polar coordinates on the z-plane: w 2 = z — r exp i<p: 

r = 1 + 2 cos </>. (150) 

It is a particular case of the trisectrix. The trivial solution holds outside the Pascal's limagon while the non-trivial 
inside. For the trivial solution the Green's function is G = Qmm = l/( z — 1) an d thus p(x,y) = 0. The non-trivial 
solution can be found by inverting Qm = {Z — TZm)^ 1 for JZm (|145[) . The Green's function is given by the upper left 
element of Gm 

r = c ^-l + fg n _ 

(^-l + |C)(^-l + fC)+C(f +2 + |) [b) 

with C being a solution of (|148[) . and again agrees with [2l|. One can write down the solution in polar coordinates 
on the z-plane: z — re 1 ^ as 



G = G X - iG y = + _ &±^*± (15 2) 



where C and D are real non-negative functions 

C = i (-1- 2r + a/1 + 8r(l +cos0)) (153) 

and 

D = ((r + C) cos (j) - l) 2 + (r + C) 2 sin 2 + 2(7(1 + cos 0) (154) 

The first one corresponds to (|148[) and the second one to the denominator in (|151l) . One can explicitly see that C is 
positive inside the Pascal limagon r < 1 + 2 cos 0. Using the Gauss law we find the eigenvalue density 

ldG 1 (8G X ^ dG y \ 1 
7T az 27r \ aa; ay / 27r 

The imaginary part of dgG is proportional to the rotation rotG = d x G y — d y G x that vanishes by construction. This 
fact can be used as a test of correctness of calculations. The density calculated from this formula is shown in figure 
ITD1 Finally, we perform some numerical checks. We generate numerically matrices M — AB = (1 + X^)(l + Xb) of 
dimensions 100 x 100 and compare obtained eigenvalue histograms with the exact solution for infinite dimensions. In 
figure [TT] we show a scattered plot of eigenvalues and the histogram of real eigenvalues compared to the section of the 
analytic solution along the real axis. The results show a good agreement between numerical data and the analytic 
result. The small remaining deviations can be attributed to finite size effects. 
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FIG. 10. The eigenvalue density of the product of two shifted Ginibre-Girko matrices. It is non-zero in the region r < 1 + 2 cos <fi- 
The density is peeked around the origin. The maximum of the function is located at the origin: p{x = 0,y = 0) = | « 1.90986 
while the minimum at the point x — 3, y = 0: p(x = 3, y = 0) = ||- « 0.0511569. 



2 |- 




FIG. 11. (Left) The analytical contour r = 1 + 2cos0 (|150[) and the scattered plot of eigenvalues obtained by diagonalisation 
of 100 matrices of dimensions 100 x 100. One should note that the boundary of the support is formed only by the external part 
of the Pascal's limacon which corresponds to £ [— 2n/3, 2n/3]. The remaining part of the limagon lies inside the support. 
(Right) A numerical histogram (solid line) constructed from almost real eigenvalues (whose imaginary part is less than e = 10 -2 ) 
obtained by diagonalisation of 20000 matrices of size 100 x 100 compared to the section of the analytic eigenvalue density p 
along the real axis. The deviations between the numerical histogram and the theoretical curve are caused by finite size effects. 

VII. SUMMARY 

We have introduced a natural generalization of the concept of S transform for the product of non-hermitian ensem- 
bles. This construction puts on the same footing addition and multiplication laws for hermitian and non-hermitian 
ensembles. We have also found a more general reformulation of the multiplication law which allows us to calculate 
free products of random matrices having vanishing mean, including the case when both factors in the product are 
centered. This case is especially interesting as it cannot be addressed using ordinary S transform techniques. 

Our construction relies on the insights from diagrammatic techniques, and in particular assumes the finiteness of 
the moments. We are however convinced, that these conditions are neither restrictive nor mandatory for a general 
proof, based on purely algebraic structures like e.g. amalgamation of free random variables and a careful treatment 
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of regularization of ensembles with ubounded moments. 
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